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Abstract 

We report a study of the electronic dissociation energy of the water dimer using quantum Monte 
Carlo (QMC) techniques. We have performed variational quantum Monte Carlo (VMC) and diffu- 
sion quantum Monte Carlo (DMC) calculations of the electronic ground state of the water monomer 
and dimer using all-electron and pseudopotential approaches. We have used Slater-Jastrow trial 
wave functions with B3LYP-like single-particle orbitals, into which we have incorporated backflow 
correlations. When backflow correlations are introduced, the total energy of the water monomer 
decreases by about 4-5 mHa, yielding a DMC energy of —76.42830(5) Ha, which is only 10 mHa 
above the experimental value. In our pseudopotential DMC calculations, we have compared the 
total energies of the water monomer and dimer obtained using the locality approximation with 
those from the variational scheme recently proposed by Casula [Phys. Rev. B 74, 161102(R) 
(2006)]. The time step errors in the Casula scheme are larger and the extrapolation of the energy 
to zero time step always lies above the result obtained with the locality approximation. However, 
the errors cancel when energy differences are taken, yielding electronic dissociation energies within 
error bars of each other. The dissociation energies obtained in our various all-electron and pseu- 
dopotential calculations range between 5.03(7) and 5.47(9) kcal/mol and are in good agreement 
with experiment. Our calculations give monomer dipole moments which range between 1.897(2) 
and 1.909(4) Debye and dimer dipole moments which range between 2.628(6) and 2.672(5) Debye. 
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I. INTRODUCTION 



Water, as the main agent of all aqueous phenomena and an important component of 
the vast majority of all chemical and biological processes, has been the subject of many 
experimental and theoretical studies. The characteristic physical and chemical properties 
of water stem from its strong polar hydrogen bonds. In spite of the apparent simplicity of 
hydrogen bonding - a hydrogen atom bonded to an electronegative group and a lone electron 
pair on another system - understanding hydrogen bonding has proved very difficult. 

One would hope to gain some insight by first investigating the water monomer and 
then progressing to the water dimer and larger water clusters, systematically including 
the many-body behaviour as the number of molecules increases. The study of the structure 
and energetics of assemblies of water molecules poses severe challenges to computational 
methods because they are held together by hydrogen bonds with binding energies of only a 
few kcal/mol. 

The hydrogen bond in the water dimer has been the subject of many electronic structure 
studies ll, 2, since it represents the prototype of all hydrogen-bonded sys- 

tems. In recent years the water dimer has been studied using high-level quantum chemistry 
techniques, such as second-order M0ller-Plesset perturbation theory (MP2) [3] and coupled- 
cluster CCSD(T) methods j^. These methods can treat electron correlation effects quite 
accurately, but they are very expensive for larger water clusters, as the costs of MP2 and 
CCSD(T) calculations scale as A^^ and A^^, respectively, where is the number of elec- 
trons. The accuracy of the results depends significantly on the quality of the basis set, and 
corrections are normally applied for the effects of basis incompleteness. Density-functional 
theory (DFT) shows a more favourable scaling with system size, allowing the use of larger 
basis sets resulting in smaller basis set errors. However, the calculated energies of the water 
monomer and dimer depend significantly on the exchange-correlation functional used. 

Quantum Monte Carlo (QMC) methods js] represent an alternative and promising way 
of treating electron correlation. The cost of calculating energies with QMC scales roughly 
as A^^, allowing accurate calculations for large systems. Moreover, QMC algorithms are 
intrinsically parallel and Monte Carlo codes are easily adapted to parallel computers. As 
the availability of powerful computers has increased, QMC has become a very attractive and 
effective method for probing the electronic structure of molecules and solids 10, 111! ]. 
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In this paper we report variational and diffusion Monte Carlo (VMC and DMC) calcu- 
lations of the water monomer and dimer, using all-electron (AE) and pseudopotential (PP) 
approaches. We have used Slater-Jastrow (SJ) and Slater-Jastrow-Backflow (SJB) wave 
functions, finding that the introduction of backflow correlations retrieves a substantial addi- 
tional fraction of the correlation energy. Our calculations yield dissociation energies in very 
good agreement with experiment. We believe that our calculations form the most accurate 
QMC study of the water monomer and dimer performed to date. 

The rest of the paper is organized as follows. In Section [TTl we give short reviews of the 
VMC and DMC methods, the SJ and SJB wave functions, and the methods used for dealing 
with non-local PPs within DMC. In Section IIIII we discuss the single-particle orbitals used 
and some details of the calculations. In Section [IV] we present and discuss our results. The 
first part of Section HVl describes results for the water monomer, the second part gives results 
for the water dimer and its dissociation energy into two water molecules, and we compare 
results from the two different schemes used to calculate the non-local PP energy. In the 
third part, we report our results for the dipole moment of the monomer and dimer. Our 
conclusions are summarized in Section |Vl 

II. VMC AND DMC METHODS 

In VMC the ground-state energy is estimated as the expectation value of the Hamilto- 
nian with an approximate trial wave function, the integrals being evaluated by importance- 
sampled Monte Carlo integration. The trial wave functions contain variable parameters, 
whose values are obtained from an optimisation procedure formulated within VMC. There 
are no restrictions on the form of trial wave functions that can be used, and VMC does not 
suffer from a fermion sign problem. However, the choice of the approximate trial wave func- 
tion is very important as it directly determines the accuracy and statistical efficiency of the 
calculation. Due to the difficulty of preparing trial wave functions of equivalent accuracy for 
different systems, the bias in the energy difference between systems is normally significant. 
We have used VMC methods mainly to optimize parameters in the trial wave functions, and 
our most accurate calculations have been performed with the DMC method. 

In DMC the ground-state component of a trial wave function is projected out by evolving 
an ensemble of electronic configurations using the imaginary-time Schrodinger equation. 
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The fermionic symmetry is maintained by the fixed-node (FN) approximation [1^ , in which 
the nodal surface of the DMC wave function is constrained to equal that of the trial wave 
function. The FN DMC energy is higher than the exact ground-state energy, and becomes 
equal to it when the fixed nodal surface is exact. The dependence of the DMC energy on 
the trial wave function is smaller than in VMC, but in practice it is often significant. It is 
therefore desirable to be able to improve the nodal surfaces in order to reduce the impact 
of the FN approximation. 

The standard SJ wave function can be written as 

vl/sj(R) = e-^(^)vl/g(R), (1) 

where R is a 3A^-dimensional vector denoting the position of each electron. The nodes 
of \I^sj(R-) are defined by the Slater part of the wave function, \E's(R), which takes the form 
\l/g = D^, where is a Slater determinant of single particle orbitals of spin a. 

The Jastrow correlation factor, e"^^^\ contains electron-electron, electron-nucleus and 
electron-electron-nucleus terms, as described in Ref. |l3|. The Jastrow factor keeps electrons 



away from each other and greatly improves wave functions in general, but it does not modify 
the nodal surface of the wave function. 

One way of reducing the FN error is to alter the nodes of the wave function by introducing 



backfiow correlations 



14| , replacing the coordinates R in the Slater part of the wave function 



by the collective coordinates X, so that the SJB wave function reads 

vl/sjB(R) = e-^(^)vI;g(X). (2) 

The new coordinates for each electron are given by 

Xi = ri + e.(R), (3) 

being the backfiow displacement of particle i, which depends on the position of every 
electron in the system. Details of the specific form of the backfiow function used can be 



found in Ref. 
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A. Pseudopotential DMC calculations 



We use non-local (angular-momentum dependent) PPs. Unfortunately non-local poten- 
tials cannot be used directly in DMC calculations because of the local character of the 

4 



algorithm. Indeed, as we have already mentioned, the starting point of the DMC algo- 
rithm is the imaginary time Schrodinger equation. If the FN Hamiltonian, H, contains a 
PP with a non-local part, V^i, and is a trial energy, then the propagator for imagi- 
nary time r is < R|exp[— r(H — Et)]|R' >. This propagator contains terms of the form 
< R|exp[— rVni] |R' >, which are not guaranteed to be positive for arbitrary R and R', and 
therefore cannot be interpreted as a probability density. 

One way of avoiding the possible change of sign is to use the PP locality approximation 
(PLA) [l^. In this approximation the FN Hamiltonian is replaced by an effective Hamilto- 
nian in which the non-local part operates on the trial wave function This technique is 
analogous to the FN approximation. In this way the non-local part of the Hamiltonian is 
replaced by a local term and the PLA effective Hamiltonian reads: 

/ dR' < R'iKilR > *t(RO 



PLA 
R,R 



K + Vioc(R 



(4) 



*t(R) 

where K is the kinetic energy operator and Vioc is the local part of the PP. In this approx- 
imation, all of the matrix elements of the localized potential enter the branching factor of 
the DMC algorithm. A disadvantage of the PLA is that the energy may be higher or lower 
than the exact value. 

Casula has recently introduced a different scheme 16|] for treating non-local PPs within 
DMC. In the Casula scheme (CS), the FN Hamiltonian is substituted by an effective Hamil- 
tonian 



R' 



^R',R=<R'IKi|R> 





rrCS 
-"R',R 



if Vn'n < 

if ^R'R > 



(5) 



where 



and 



R'R 



^t(R') 
*t(R) 



< R'iKilR > 



1 



(6) 



(7) 



From Eqs. ([5]) and ([6]) we see that in the CS only the positive matrix elements of the non- 
local potential are localized, and they are absorbed into the branching factor. On the other 
hand, unlike in the PLA, the negative matrix elements of Vai contribute to the drift-diffusion 
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of the walkers. In practice, the standard DMC algorithm needs few changes: once a drift- 
diffusion move is proposed and accepted or rejected, and after the walker has been weighted, 
an extra displacement of the walker is performed according to a transition probability that 
depends only on the negative matrix elements of the non-local part of the PP. This method 
has two advantages over the PLA. First, the ground state energy of the CS Hamiltonian 
is an upper bound [l7| on the ground state energy of the true FN Hamiltonian. Second, 
whenever the negative elements are large, the new displacement pushes the walkers away 
from the attractive regions of the localized potential, which avoids "population explosions" 
in which walkers are trapped in the attractive regions with large branching factors. We 
found occasional population explosions in our PLA calculations, which we dealt with by 
restarting the calculations at an earlier move with a new random seed. We found no such 
instabilities when using the CS. 



III. DETAILS OF THE CALCULATIONS 



All the calculations have been performed using the casino code [l^. The single-particle 
orbitals forming the Slater determinants have been obtained from DFT calculations using 
the CRYSTAL98 code [l9| with the Roos augmented double zeta ANO (Roos aug-DZ-ANO) 
Gaussian basis set j20| used in its decontracted form. Extensive tests of the basis sets 



listed at http://www.emsl.pnl.gov/forms/basisform.html showed that the Roos basis set 



gives very good DMC energies. The dependence of the DMC energies on the orbitals is 



generally quite weak. However, previous calculations for neon, neon radicals 



211 ] and water 



22( 1 have shown that B3LYP orbitals give slightly lower DMC energies than Hartree-Fock 



(HF) orbitals, suggesting that the former can improve the nodal surface of the wave function. 
Based on this, we have investigated whether one can obtain further improvements in the 
DMC energies by optimising the form of the exchange-correlation (XC) functional used to 
generate the single-particle orbitals. Becke's three parameter B3LYP XC functional 
may be written as 19|: 



23| 



1 - A){E^^^ + B E™) + A E^"' + (1 - C) E^'''' + CE, 



nBecke 



HF 



r-VWN 



r-LYP 



where A determines the amount of Fock exchange, and B and C determine the amounts of 
non-local exchange and correlation, respectively. In Eq. (j8]), the Vosko-Wilk-Nusair corre- 
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lation potential (VWN) is used to extract the local part of the LYP correlation potential. 
The usual B3LYP functional corresponds to A = 0.2, B = 0.9 and C = 0.81. In principle we 
could try and minimize the DMC energy with respect to A, B and C, but this would be very 
expensive. Setting ^4 = 1 in Eq. ([S]), we have obtained a DMC energy for the water molecule 
of -76.4218(1) Ha [44j, which is very close to the value of -76.42205(8) Ha obtained with 
HF orbitals. This suggested that the contribution from the correlation part (terms involv- 
ing C) does not significantly affect the DMC energy. We therefore decided to search for the 
amount of Fock exchange {A) which minimizes the DMC energy, while keeping B and C 
equal to their standard B3LYP values. The results are shown in Fig. [TJ Note that the DMC 
energy (squares) is much less sensitive to the changes in the XC functional than the DFT 
energy (circles). The dashed line represents the best quadratic fit to the DMC energies for 
A <0.6. The lowest DMC energy is obtained when 25% Fock exchange is used in the above 
functional. We therefore used A = 0.25 to generate the single-particle orbitals for the Slater 
parts of all of the trial wave functions used in this work. In what follows, we refer to these 
as B3LYP-like orbitals. 



In our AE calculations, the single-particle mo 



corrected using the scheme described in Ref. 



24 



ecular orbitals close to the nuclei have been 
so that each orbital obeys the appropriate 



cusp condition 25| at each nucleus, ensuring that the local energy is finite at each nucleus. 



In our PP calculations, we have used the two schemes described in Section I 



3 with PPs 



27|. These 



generated within Hartree-Fock theory, including scalar relativistic effects 126 
PPs have been found to work very well in conjuction with QMC methods j21[. One of the 
advantages of using PPs is that they avoid the short-range variations in the wave function 
near the nuclei, and hence we can use larger time steps than is possible in AE calculations. 
We have performed calculations for different time steps and then extrapolated to the limit of 
zero time step by fitting to quadratic functions. Population-control biases have been checked 
by running with different population sizes, and no such bias was detected. 

The parameters in the Jastrow and backfiow functions were obtained by minimizing the 
variance of the local energy 28|, l29[, except for the AE calculations with SJB wave functions. 
In this case, the Jastrow and backfiow parameters were obtained by minimizing the mean 
absolute deviation of the set of the local energies from the median. 
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IV. RESULTS AND DISCUSSION 



The calculations for the water monomer were carried out in the experimental equilibrium 
geometry with ron = toh' = 0.9572 A and L = 104.52°. This geometry was used in previously 
reported QMC calculations [22^, yOl • 

n 

For the water dimer we have chosen the geometry derived by Klopper et al. [5| from 
CCSD(T) calculations, extrapolating to the basis set limit using MP2-R12 results. The op- 
timized water dimer geometry (see Fig. [2]) yields a final equilibrium oxygen-oxygen distance 
of 2.912 A. When the dimer is formed, some structural deformation occurs in each of the 
monomers: a slight change in the angle between the inter-oxygen vector and the plane of 
the proton-accepting monomer occurs and the bond angle for the proton acceptor widens 
by about 0.5°. We have neglected these small changes, assuming that their effect on the 
intermolecular binding energy is small. The electronic dissociation energy in the equilibrium 
geometry, Dg, with respect to dissociation into two isolated monomers, is defined as 

Dg = Do + ZPE(dimer) - 2 ZPE(monomer), (9) 

where Do is the experimental dissociation energy and ZPE denotes the zero point energy of 
the nuclei. 

A. Water monomer 

Our AE SJ and SJB QMC results with B3LYP-like orbitals for the water monomer are 
given in Table [U together with energies calculated using other ab initio methods. The 
expansion order of the polynomials and the spin dependencies used in the Jastrow factor as 



described in Ref. 



13l were the following: N,i=N^=6, Nj^''^'^=3, Su=l and S;^=S/=0. This led 



to 128 parameters in the Jastrow factor. When backfiow correlations were incorporated the 
total number of parameters to be optimized increased to 338, corresponding to the following 



expansion order of the polynomials in the backfiow term (Ref. Il4l ): N^=N^=6, N^^'^'^=3, 
S,;=S0=1 and S^=0. We have mentioned in Section [III that the VMC energy depends strongly 
on the trial wave function. This is illustrated by Table [H when backfiow correlations are 
included, the VMC energy for the water molecule [—76.4034(2) Ha] decreases by nearly 
26 mHa. Our DMC energy extrapolated to zero time step with the SJ wave function (DMC- 
SJ), —76.42376(5) Ha, is nearly 1 mHa lower than that recently reported by Benedek et 
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al. 



22| using true B3LYP orbitals and lies about 14 mHa above the experimental result of 



—76.438 Ha. When backflow correlations are introduced, the DMC energy decreases by a 
further 4.5 mHa, which is 25% of the correlation energy missing with the SJ wave function, 
yielding a DMC energy of —76.42830(5) Ha. Note that, as mentioned in Section [TTl the 
difference between the energies obtained with the SJ and SJB trial wave functions is smaller 
at the DMC level than at the VMC level. 

Liichow and Fink 3l| performed DMC calculations using a pair natural orbital 



configuration-interaction (PNO-CI) trial wave function consisting of 300 determinants and 
a Jastrow factor, obtaining an energy of —76.429(1) Ha, which is within error bars of our 
single-determinant SJB DMC result. For completeness. Table |T] also gives the DMC energy 



obtained by Casula et al. 



32j using a trial wave function consisting of a product of an 



antisymmetrized geminal power (AGP) wave function and a Jastrow factor, which is nearly 

so gives the DMC energies calculated with 



22 



As was pointed out in Section 



7 mHa above our best DMC result. Table m ali 
SJ wave functions and HF nodes from Refs. l2l| and 
both of these energies are slightly higher than the DMC results obtained using B3LYP and 
B3LYP-like orbitals, which suggests that the hybrid DFT functional yields better nodes than 
the HF orbitals. The small difference between the two DMC-HF values given in Tabled arises 
from the use of different basis sets. For comparison, we have also included energies from a 
CCSD(T)-R12 calculation 33| and a recently reported Correlation Energy Extrapolation by 
Intrinsic Scaling (CEEIS) calculation 3J] which, although not variational, give results very 
close to the "exact" one. 



B. Water dimer and electronic dissociation energy 

As illustrated in Fig. [21 the two hydrogen atoms of the acceptor monomer (Ha) are not 
equivalent to the donor hydrogen (Hd) and free hydrogen (Hf) atoms in the donor subunit, 
and the two oxygen atoms (Oa and Od) are inequivalent. However, we have obtained slightly 
lower DMC energies by using Jastrow and backflow functions which treat all four hydrogen 
atoms and both oxygen atoms as equivalent. Treating atoms as equivalent in this way reduces 
the number of variable parameters, and presumably this results in the better performance of 
the stochastic optimization procedure. When the H and O atoms are treated as equivalent, 
the SJ trial wave functions for the dimer contains 128 parameters (for both the AE and PP 
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calculations) while the SJB calculations have been performed with 390 parameters for the 
AE calculations and 286 parameters for the PP calculations. These number of parameters 
correspond to N„=N^=6, Nj^'^'^=3, S„=S/=1 and S^=0 in both the AE and PP Jastrow 
factors. In the backflow term, N^=N^=6, S^=S^=1, S^=0 in both AE and PP calculations, 
and N^''''^"=3 and N^'''""=2 in the AE and PP cases, respectively. 

Figure [3] shows the AE DMC energies of two water monomers and the water dimer as 
a function of time step. The zero-time-step extrapolations of the fitted curves are shown 
by filled triangles of the same color as the original symbols. The DMC-SJ and DMC-SJB 
values shown in Table [T] correspond to the black and blue triangles in this figure, i.e., the 
extrapolated values from the curves labelled 2xH20(SJ) and 2xH20(SJB), respectively. 
As we mentioned previously for the water monomer, the DMC energy decreases by about 
4.5 mHa when backflow correlations are introduced. Figure [3] shows that the reductions in 
energy from including backflow correlations are essentially independent of the time step for 
both the monomer and dimer. 

Figure m shows PP SJ DMC energies for two water monomers and the dimer as a function 
of the time step. The CS and PLA curves cross at a time step of about 0.02 au. The 
extrapolated SJ DMC energy of two water monomers obtained within the CS scheme is less 
than 4 mHa above the PLA value. It is not surprising that the CS energies are above the 
PLA ones, because the CS scheme gives an upper bound on the true energy, while the PLA 
scheme does not. The energy of the water monomer as a function of the time step in both 
the PLA and CS schemes can be found in Table [Tll 

In Fig. [5] we have summarized the total energies of the PP water monomer and dimer 
energies with SJ and SJB wave functions within the PLA. The CS yields a similar figure (not 
shown) with curves parallel to the red squares of Fig. HI As in Fig. [3l including backflow 
correlation in the monomer and dimer results in an essentially rigid shift of the corresponding 
SJ curve to lower energies. Furthermore, because the black and red curve corresponding to 
the SJ calculations, and the blue and green curves, corresponding to the SJB calculations, 
have nearly the same shape, the time step errors largely cancel when computing the De of 
the water dimer. This cancellation is illustrated in Fig. [6l where we have plotted D^, as a 
function of the time step for the four combinations of SJ/SJB wave functions within the 
AE/PP approaches from Figs. [3] and [5l 

The top panel of Fig. [6] shows the AE values of D^. The magenta diamonds (orange 
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triangles) show the results with SJ (SJB) wave functions as a function of time step. We 
fitted the energies to quadratic functions and extrapolated to zero time step, obtaining the 
values of D^, given in Table lllli 

The high degree of cancellation of the time step errors shown in Fig. E] is remarkable. Our 
AE-SJ data show smaller time step errors than those reported in Ref. |22| with HF single- 
particle orbitals. When backfiow correlations are introduced, the cancellation of time step 
errors at larger time steps is not as good. This can also be observed in Fig. [3] where the 
green line bends more than the blue one at larger time steps. Nevertheless, the extrapolated 
AE values of for the SJ wave function [5.16(10) kcal/mol] and the SJB wave function 
[5.35(10) kcal/mol] are within error bars of each other. For the PP calculations we have seen 
that the PLA and CS schemes give significantly different time step errors. However, the time 
step errors are very similar for the monomer and dimer (see Fig. Hj). Therefore the values 
of Df, computed within the PLA and CS schemes are very similar at each time step. We 
find De = 5.03(7) kcal/mol [5.07(7) kcal/mol] for SJ wave functions within the PLA [CS], 
and De = 5.47(9) kcal/mol [5.38(8) kcal/mol] when backfiow correlations are introduced. 
We have also found that, in both our AE and PP calculations, backfiow correlations tend 
to increase the value of slightly. 

Our DMC results range between 5.03(7) and 5.47(9) kcal/mol and are similar to the DMC 



values of Df. reported in Ref. |22| (with HF and true B3LYP orbitals) and the CCSD(T) result 



of Klopper et al. [5l] of 5.02(5) kcal/mol. 

For many years the experimental value of has been taken to be 5.44±0.7 kcal/mol. 
This value was derived from measurements of the enthalpy of dimerization of water corrected 
with a theoretical estimate of the ZPE 36|] . Mas et al. 37| have recently calculated a more 



accurate value of the ZPE and have estimated to be 5.00±0.7 kcal/mol. Leforestier et 



al. [38| have reported a value of Df, of 5.14 kcal/mol, correcting for the ZPE by determining 
potential energy surfaces for the monomer and dimer via direct inversion of spectroscopic 
data. 

C. QMC dipole moment of the water monomer and dimer 

We have calculated the mean dipole moment, /i, of the water monomer and dimer. The 
dipole moment operator does not commute with the Hamiltonian and therefore the error in 
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the DMC mixed estimate (/xdmc) is linear in the error in the wave function. This error can 
be reduced to a quadratic error by using the so called extrapolated estimator, 2/iDMC — /Uvmc, 
where /ivMC is the VMC estimate of the dipole moment. Both VMC and DMC estimates 
must, of course, be calculated using the same trial wave function. 

Figure [7] shows our extrapolated estimator results for the dipole moment of a water 
molecule as a function of time step. Red squares (black triangles) show AE SJ (SJB) results 
while green circles (blue diamonds) represent PP SJ (SJB) results calculated with the CS. 
Dashed lines of the same color as the original symbols are linear fits used to find the zero- 
time-step dipole moment. Our results are summarized in Table HVl together with the results 
of other ab initio calculations. Our results for the monomer are slightly above (0.05 Debye) 
the experimental value, as are the other theoretical values. For comparison, the inset of 
Fig. [7] shows the DMC mixed estimates of the dipole moment for each of the four cases 
investigated. As we have mentioned above, the DMC estimate has a larger dependency on 
the trial wave function and this can be observed in the larger differences between the data 
sets displayed in the inset. However, when the extrapolated estimator is considered the error 
is reduced and the results come into closer agreement. 

We have arrived at similar conclusions for the dipole moment of the water dimer (figure 
not shown). The extrapolated estimator values of the dipole moment of the water dimer, 
extrapolated to zero time step, are summarized in Table |Vl Our AE results are lower than 
the experimental value 39[, but within 0.02 Debye of it, and our PP results are higher than 
the experimental value, but within 0.03 Debye of it. 



V. CONCLUSIONS 



We have performed all-electron and pseudopotential QMC calculations of the total energy 
of the water monomer and dimer using SJ and SJB wave functions with B3LYP-like single- 
particle orbitals. The introduction of backfiow correlations recovers a significant amount 
of correlation energy and we have obtained a final DMC energy (extrapolated to zero time 
step) for the water molecule of —76.42830(5) Ha. We have performed a careful study of the 
time step errors and have shown that they largely cancel when energy differences, such as 
the electronic dissociation energy of the water dimer, are considered. We have compared two 
different schemes for treating non-local PPs within DMC: the pseudopotential localization 
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approximation and a scheme recently proposed by Casula. Although the total energies 
from the two schemes differ significantly, these differences are essentially the same in the 
monomer and dimer, and therefore they cancel when evaluating D^. Our extrapolated AE 
and PP results for the D^. show good agreement with each other when both SJ and SJB 
wave functions are used, and they are all within the error bars of the values determined 
from experiment. The total spread of our AE and PP values of is only 0.44 kcal/mol. 
We have also reported QMC estimates of the dipole moments of the water monomer and 
dimer. Our AE and PP extrapolated DMC results yield a dipole moment of the monomer 
within 0.05 Debye of the experimental value, while for the dimer the AE results are within 
0.02 Debye of experiment and the PP results are within 0.03 Debye of experiment. 
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Figure captions 



Fig. 1: DMC energy (squares) as a function of the percentage of Fock excliange in the 
exchange- correlation B3LYP-Iike functional [A in Eq. ([H])] used to obtain the single-particle 
orbitals for the SJ wave functions. The same Jastrow factor of 128 parameters was used in all 
the calculations. The dashed line shows a quadratic fit to the points around the minimum. 
The dotted line represents the DMC energy with HF orbitals. Inset: DFT energy as a 
function of the percentage of Fock exchange. 



Fig. 2: (Color online) The equilibrium structure of the water dimer from Klopper et al. 
jsi The "a" stands for acceptor, "d" for donor and "f" for free. 



Fig. 3: (Color online) All-electron total energies of two water monomers and the wa- 
ter dimer as a function of time step, calculated with SJ (circles) and SJB (squares) wave 
functions. Dashed and dotted lines show quadratic fits to the monomer and dimer energies, 
respectively. The filled triangles represent the extrapolation to zero time step. The statisti- 
cal error bars are of order 8 x 10^^ Ha. 



Fig. 4: (Color online) SJ PP energies for two water monomers and for the water dimer, 
calculated with the PLA (black circles) and the CS (red squares). Dashed (dotted) lines are 
fits for the monomer (dimer). Triangles denote the energies extrapolated to zero time step. 
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Fig. 5: (Color online) As in Fig. [3], but using PPs within the PLA. 



Fig. 6: (Color online) Dissociation energy of the water dimer as a function of time step. 
The curves represent time step extrapolations and filled symbols give the zero time step D^. 
Top panel: AE calculations. Bottom panel: PP calculations. For comparison, the AE values 
extrapolated to zero time step are also shown (with the same colors and symbols as in the 
top panel). 



Fig. 7: (Color online) Extrapolated estimates of the dipole moment of the water molecule 
as a function of the time step. Filled symbols correspond to zero time step values obtained 
from linear fits (dashed lines). Inset: DMC mixed estimator. 
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Method 


Total energy (Ha) 


Basis set 


Ref. 


DFT 


-76.4361 


Roos Aug-DZ-ANO 


This work 


VMC-SJ 


-76.3773(2) 


Roos Aug-DZ-ANO 


This work 


VMC-SJB 


-76.4034(2) 


Roos Aug-DZ-ANO 


This work 


DMC-SJ 


-76.42376(5) 


Roos Aug-DZ-ANO 


This work 


DMC-SJB 


-76.42830(5) 


Roos Aug-DZ-ANO 


This work 


DMC-HF 


-76.42102(4) 


6-311+-FG(2d,2p) 


21 


DMC-HF 


-76.4219(1) 


Roos Aug-DZ-ANO 


22 


DMC-B3LYP 


-76.4230(1) 


Roos Aug-DZ-ANO 


22 


DMC-PNO-CI 


-76.429(1) 


STO 


31 


DMC-AGP 


-76.4175(4) 


STO-DZ 


32 


CCSD(T)-R12 


-76.4373 


Dunning cc-pV6Z 


33 


CEEIS 


-76.4390(4) 


CBS limit 


34 


Exact 


-76.438 




35 



TABLE I: Total energy of the water monomer calculated with different ab initio methods, and the 
exact non-relativistic static point nuclei value. CBS limit stands for extrapolation to the complete 
basis set limit. 
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Energy (Ha) 
Time step (au) PLA CS 

0.000 -17.20843(4) -17.20678(4) 

0.005 -17.20821(3) -17.20698(3) 

0.010 -17.20791(3) -17.20714(3) 

0.020 -17.20762(3) -17.20763(3) 

0.030 -17.20757(3) -17.20832(2) 

0.050 -17.20783(2) -17.20994(2) 



TABLE IL Comparison of PP SJ DMC energies obtained within the PLA and CS schemes. 
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Method 


Df, (kcal/mol) 


Ref. 


AE-SJ 


5.16(10) 


This work 


AE-SJB 


5.35(10) 


This work 


PP-SJ (PLA) 


5.03(7) 


This work 


PP-SJ (CS) 


5.07(7) 


This work 


PP-SJB (PLA) 


5.47(9) 


This work 


PP-SJB (CS) 


5.38(8) 


This work 


AE-HF-SJ 


5.02(18) 


22 


AE-B3LYP-SJ 


5.21(18) 


22 


CCSD(T) 


5.02(5) 




Experiment + Theory 


5.44(70) 


36 


Experiment + Theory 


5.00(70) 


37 


Experiment + Theory 


5.14 


38 



TABLE III: Top: DMC dissociation energies (kcal/mol) obtained in this work with SJ and SJB wave 
functions using B3LYP-like orbitals within the AE and PP (PLA and CS) approaches. Bottom: 
Results from other ah initio calculations and from "Experiment + Theory". No error bar was 



reported in Ref. 



38 
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Method (Debye) Ref. 

AE-SJ 1.909(4) This work 

AE-SJB 1.900(3) This work 

PP-SJ (CS) 1.897(2) This work 

PP-SJB (CS) 1.900(3) This work 

aug-cc-VDZ MP2 1.868 40 

aug-cc-pVQZ MP2 and MP4 1.87 41 
aug-cc-pVDZ MP2-R12 1.883 42 
Experiment 1.855 40 

TABLE IV: Top: Dipole moment of the water monomer, evaluated using the extrapolated estima- 
tor, and also extrapolated to zero time step. Bottom: Other ab initio results and the experimental 
value. 
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Method 


(Debye) 


Ref. 


AE-SJ 


2.628(6) 


This work 


AE-SJB 


2.635(7) 


This work 


PP-SJ (CS) 


2.670(5) 


This work 


PP-SJB (CS) 


2.672(5) 


This work 


aug-cc-VDZ MP2 


2.683 


40 


aug-cc-pVDZ B3LYP 2.635 


43 


Experiment 


2.643 


39 



TABLE V: As in Table IVl but for the water dimer. 
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